#include <assert.h>
#include <stdio.h>

#include "../src/mstep.c"

typedef struct _MstepTestDatum MstepTestDatum ;
struct _MstepTestDatum 
{ 
  struct
    {
      float     kappa;
      float     logalpha;
      float     tau[(5-1)*4];
      float     loglambda[4];
      float     pz[5];
      Rating    ratings[8];
      float     logalphastddev;
    }                   inputs;

  struct
    {
      float     newq;
      float     logalpha;
    }                   desired_output;
};

static void
test_mstep (void)
{
  MstepTestDatum data[] = 
    {{{0.0582138,  0.46808, {0.651832, 0.47304, 0.654126, 0.0216862, 0.820111,  0.101254, 0.340756, 0.0672765, 0.0541735, 0.427385, 0.0324973,  0.643789, 0.761232, 0.215306, 0.328076, 0.756206}, {0.0971072,  0.863713, 0.31651, 0.80352}, {0.0292678, 0.00453248, 0.258297,  0.335441,  0.377435}, {{0, 1}, {0, 0}, {0, 0}, {1, 2}, {1, 3}, {2, 0}, {0,  4}, {0, 4}},  0.531492}, {-17.519295827205208997, -1.1503294729550958020}},  {{0.604171,  0.313754, {0.557324, 0.430238, 0.263415, 0.246478, 0.503151,  0.0028535, 0.230917, 0.602689, 0.741919, 0.787548, 0.902841,  0.846482, 0.644812, 0.923835, 0.586331, 0.0429621}, {0.615544,  0.919303, 0.328034, 0.707522}, {0.238108, 0.38781, 0.723863,  0.393767,  0.680784}, {{3, 1}, {0, 0}, {0, 1}, {2, 0}, {2, 3}, {1, 0}, {1,  4}, {2, 2}},  0.957572}, {-37.179865943551341673, -2.2983501638432758765}},  {{0.460448,  0.147289, {0.177633, 0.954719, 0.229531, 0.5446, 0.435714,  0.167171, 0.32669, 0.698118, 0.790902, 0.243336, 0.740359,  0.655156, 0.175358, 0.324034, 0.412325, 0.947634}, {0.93725,  0.936223, 0.688462, 0.553867}, {0.256466, 0.978651, 0.228013,  0.406578,  0.0788331}, {{2, 2}, {1, 4}, {1, 2}, {3, 3}, {0, 4}, {3, 4}, {0,  3}, {2, 2}}, 0.0239323}, {-149.39131693785934473,  0.38172677157502423460}}, {{0.998483,  0.861978, {0.643119, 0.856761, 0.671793, 0.16386, 0.852217,  0.613425, 0.931434, 0.508704, 0.676858, 0.289392, 0.519109,  0.561069, 0.739608, 0.353168, 0.830647, 0.00720189}, {0.483142,  0.374517, 0.602634, 0.600624}, {0.404309, 0.350585, 0.604151,  0.738646,  0.76119}, {{2, 0}, {0, 1}, {0, 1}, {0, 2}, {3, 4}, {3, 3}, {3,  2}, {3, 1}},  0.493824}, {-51.370989476775160557, -1.3306325324662280198}},  {{0.932358,  0.574786, {0.908973, 0.880399, 0.000924218, 0.0660824, 0.232115,  0.591007, 0.481815, 0.505013, 0.492507, 0.237838, 0.651168,  0.497811, 0.00936459, 0.863321, 0.0485336, 0.897188}, {0.605055,  0.512736, 0.444382, 0.158542}, {0.843865, 0.0189124, 0.512024,  0.583756,  0.934892}, {{3, 4}, {3, 0}, {0, 1}, {0, 3}, {0, 3}, {3, 1}, {3,  1}, {0, 1}},  0.138514}, {-106.13496232850602594, -0.14227277311979856798}},  {{0.5111,  0.517673, {0.702777, 0.547507, 0.0292848, 0.0126602, 0.21027,  0.309669, 0.378117, 0.514849, 0.200905, 0.446347, 0.329584,  0.617661, 0.59585, 0.933611, 0.885201, 0.45912}, {0.751985,  0.914699, 0.373177, 0.875364}, {0.817093, 0.776185, 0.862077,  0.357691,  0.114316}, {{0, 0}, {1, 3}, {1, 3}, {3, 2}, {2, 1}, {1, 1}, {0,  1}, {1, 0}},  0.228678}, {-64.654880571983433656, -0.77744760889869077500}},  {{0.832793,  0.34503, {0.904046, 0.91901, 0.454675, 0.830182, 0.703141,  0.472662, 0.125092, 0.21252, 0.107291, 0.539051, 0.23989,  0.753401, 0.355306, 0.624352, 0.866713, 0.878037}, {0.538214,  0.848167, 0.00463583, 0.520346}, {0.423897, 0.619489, 0.171843,  0.175316,  0.519851}, {{1, 2}, {1, 1}, {3, 4}, {2, 1}, {2, 2}, {3, 0}, {1,  0}, {0, 4}},  0.700479}, {-34.626612076959586460, -1.9024007359987113119}},  {{0.717168,  0.345134, {0.816711, 0.227817, 0.592076, 0.132614, 0.70942,  0.688766, 0.352186, 0.379213, 0.354113, 0.0644142, 0.485473,  0.501177, 0.8159, 0.216247, 0.480837, 0.98083}, {0.392002,  0.596759, 0.308994, 0.805515}, {0.872151, 0.896279, 0.591826,  0.460381,  0.0554404}, {{2, 2}, {1, 3}, {1, 0}, {2, 2}, {0, 2}, {1, 0}, {3,  3}, {0, 1}},  0.668462}, {-42.656701901515194570, -1.5356808168417184298}},  {{0.999749,  0.327767, {0.346021, 0.979696, 0.647563, 0.948553, 0.991907,  0.915282, 0.162091, 0.447377, 0.176007, 0.699035, 0.681254,  0.466547, 0.784005, 0.102276, 0.37226, 0.661032}, {0.911854,  0.205997, 0.780435, 0.200651}, {0.856414, 0.537534, 0.780686,  0.872884,  0.510393}, {{1, 2}, {1, 0}, {3, 1}, {1, 4}, {0, 0}, {3, 1}, {3,  1}, {3, 4}},  0.557838}, {-60.661756987728827376, -1.6020699932411985382}},  {{0.133122,  0.924331, {0.518486, 0.642556, 0.971031, 0.476954, 0.342478,  0.943521, 0.289777, 0.0104076, 0.558473, 0.841245, 0.917517,  0.349376, 0.646619, 0.635249, 0.137082, 0.148725}, {0.790206,  0.0977148, 0.356397, 0.27584}, {0.279813, 0.539877, 0.223274,  0.351509,  0.761327}, {{1, 0}, {2, 2}, {0, 0}, {3, 1}, {0, 2}, {2, 4}, {3,  4}, {3, 4}},  0.897321}, {-31.984105816209243525, -2.1589694916016017479}}} ;

  for (unsigned int i = 0; i < sizeof (data) / sizeof (data[0]); ++i)
    {
      ParameterizedNormalDistribution prior_logalpha = 
        parameterized_normal_distribution (data[i].inputs.logalphastddev);

      float hq = 0;
      for (unsigned int k = 0; k < 5; ++k)
        {
          hq -= data[i].inputs.pz[k] * log (data[i].inputs.pz[k]);
        }

      float logalpha = data[i].inputs.logalpha;
      float q = mstep (data[i].inputs.kappa,
                       &logalpha,
                       5,
                       data[i].inputs.tau,
                       data[i].inputs.loglambda,
                       data[i].inputs.pz,
                       data[i].inputs.ratings,
                       8,
                       &prior_logalpha.base,
                       hq,
                       1000,
                       1e-7);

      assert (fabs (q - data[i].desired_output.newq) <=
              1e-3 * (1 + fabs (q) + fabs (data[i].desired_output.newq)) ||
              ( fprintf (stderr, "%g ?= %g\n", q, data[i].desired_output.newq),
                0 ));
      
      assert (fabs (logalpha - data[i].desired_output.logalpha) <=
              1e-3 * (1 + fabs (logalpha) + fabs (data[i].desired_output.logalpha)) ||
              ( fprintf (stderr, "(oops ... but this is ok) %g ?= %g\n", logalpha, data[i].desired_output.logalpha),
                1 ));

    }
}

int 
main (void)
{
  srand48 (69);

  test_mstep ();

  return 0;
}
